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ABSTRACT 


Various computer-oriented numerical schemes are 
available to evaluate the velocity potentials for 
inviscid fluid flow past submerged objects. The use of 
a surface distribution of sources of uniform strength 
cver panels representing the object's surface reguires 
a relatively fine grid to obtain accurate results. A 
new variation of this scheme using distributed sources 
of linearly varying strengths over triangular panels 
appears to.reduce the number of panels needed for gcod 
results, thereby reducing computer time and storage 
requirements. In the present study, this technique is 
applied tc the problem of determining fluid forces and 
dynamic response resulting from interaction between an 


arbitrarily shaped object and surface waves. 
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Ocean waves interacting with large marine structures 
preduce substantial hydrodynamic forces. These forces are 
experienced ty floating objects such as ships or moored 
platforms as well as submerged bodies near the free surface, 
such: as sukmarines at periscope depth or bottom mounted 


caissons in shallow water. 


It is desirable from an engineering standpoint to be 
able to determine the wave forces and resulting dynanic 
response for any large ocean structure under design. Until 
recently, however, such determinations were not possible in 
the general case. Application of the well-known Mcrison 
eguation to calculate wave forces is limited to cases in 
which structure dimensions are small ccmpared to wave 


length, as in the case of a small diameter pile. 


The usual analytical approach in dealing with the 
hydrodynamic interaction of large bodies with waves is based 
on a mathematical modeling of the fluid flow using potential 
flcw thecry, which assumes the fluid to te inviscid and 
irrotational. While such an assumption is unrealistic for 
steady flow rast bluff bodies (because of boundary layer 
separaticn effects), experiments indicate that it yields 
valid results for the case of unsteady flow past an cbject 
as produced by surface waves provided wave amplitude is 
small compared to the characteristic dimensions cf the 
D ject. 


Numerical schemes based on the use of the Green's 


function (cr source distribution) method which utilize 





high-speed computer capabilities have had considerable 
success in solving problems of wave/body interactions. 
Garrison, et. al. (Refs. 1-5] have developed numerical 
procedures for calculation of wave interaction with large, 
three-dimersicnal, fixed and floating bodies. These 
procedures, which represent the immersed surface by a grid 
of distributed source panels of uniform strengths, have 
proven to be powerful Ene practical applicatien as 


demonstrated in Refs. 6-9. 


The most severe limitaticn of the above mentioned 
prccedures is the rather large amount of computer (C.P.U.) 
time required to obtain an accurate solution . Thus, any 
improvements to the method which would imprcve convergence 


of the sclution would be of considerable value. 


A natural extension of the technique based on a uniform 
distribution cf sources over each grid panei is to allcw the 
source strength to vary linearly over the extent cf each 
panel. Since this is not convenient when uSing quadrilateral 
panels, it is appropriate to divide the source surface into 
a grid of trianguiar panels. The linearly varying source 
strength may then te represented in terms of the values of 


the source strength at the ccrners of the triangles. 


Webster [Ref. 10] has developed such a scheme for the 
case of uniform flow of an unbounded fluid past a 
three-dimensicnal body. The body was represented by use of 
the three-dirensional source potential of the form 1/R where 
R denotes tke distance from the source. In the present study 
this general method is extended to the problem of 
representing a body of arbitrary shape located in fluid of 


finite depth in the presence of a free surface. 





The basic approach used in solving the problem of wave 
interaction with large floating objects is to decompose it 
into compcnent parts which are more simply analyzed, andá 
then superimpese the results to yield the overall solution. 
A group cf boundary-value problems are generated in terms of 
velocity potential; one problem is solved for each of the 
six possible degrees of freedom of body moticn, and one for 
the scattering of the incident wave by the body. The 
solutions are then formulated in terms of the Green's 
function. The resulting integral equation is then converted 
to a system of linear equations which are subsequently 


sclved by ccuaputer. 





II. FORMULATION OF THE WAVE INTERACTION PROBLEM 


AÀ. RESTRICTIONS AND ASSUSPTLIONS 


All of tke usual simplifying assumptions incorporated in 
linear gravity wave theory are adopted in the present study. 
In summary, these are: 

a The flow is considered to be incompressible, 
Mmerotaticnal, and inviscid. 

" The wave amplitude to wave length ratic is 
small. 

a The pressure is uniform above the free surface, 
and the density is uniform throughout the fluid. 

" The bottom boundary is represented by a 
horizontal, impermeable plane. 

a The velocity-squared term in Bernoulli's 
equation is neglected. 


s the waves are regular and periodic. 
ii addition, it is assumed that 
s The object has a smooth, rigid, and impermeable 
surface. 


'" The wave amplitude is small compared to the 


object's dimensions. 


B. PROBLEM LESCRIPTION 


The general description of the problem is shown in Fig. 
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1, in which two ccordinate systems are defined. The inertial 
coordinate system, positioned at the mean free surface, has 
directions defined by x, y and z, while the coordinate 
system positioned on the object has directions indicated by 
x", y' and z'. The origin of the body coordinates is taken 
to be at the center of gravity of the body, and the body is 
surtmerged an arbitrary distance, d, below the mean free 
surface. Ihe incident wave is assumed to propagate in the 


direction defined by the angle X indicated in Fig. 1. 


The analysis of the flcating object located in a train 
of regular waves is decomposed into seven problems which are 
uncoupled and can be solved independently and superimposed 


to cbtain the complete solution. 


The first six problems are identical in fora, and 
consist cf determining the fluid forces đue to oscillation 
of the object in the six degrees of freedom, in an ctherwise 
still fluid. The six modes of oscillation are all considered 
to take place with angular frequency, @=2"/T, which matches 
the frequency cf the incident waves. These harmonic motions 


of the object may be expressed as 


EN 1 
X,= Re LX, e ed Ii ese = 


o o o 
in which x cae. when k=4, 5, 6, and where = denotes the 


angular displacement about the x', y! and Ze axes, 


respectively. The subscripts k=1, 2, 3 correspond to surge, 
heave, and sway modes of motion, respectively, while the 
subscripts k=4, 5, 6 correspond to the rotational modes of 
roll, yaw, and pitch, respectively, as depicted in Fig. 2. 


O 
The parameter X denotes the complex amplitude of the linear 
k 
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O 
modes of motion, while the tern 28), denotes the complex 


amplitude of the rotational degrees of freedon. 


The seventh problem to be considered is the interaction 
of waves with the object held fixed in position. This 
interaction cf an incident wave with the fixed body .results 
in "scattered" waves which radiate from the body. From the 
sclution of this problem the excitation forces and  mcments 


acting on the body are determined. 
ee EQUATICNS OF MOTION 


for a flicatine body with Linear elastic ccnstraints the 
eguation of motion for the six degrees of freedom is given 


in Ref. 1 as 


6 f P. T 
= ET aoe .. e. t -— è t i 
2 (my; eM. .) LY Ni; + (K: +K) g^ = [c;]e t (2) 
EG 
where V= ot , and 0O=2W/r denotes the frequency of the 
incident wave and resulting motion. (an PEGs 2 and all 


equations which follow, the use of tne lower case letter, i, 
as a sukscript denotes an index. Appearance of i elsewhere 
in an equaticn denotes the imaginary number i=J-1 .) The 
surface elevation of the incident wave at the origin is 


given by 


N= N cos (ct) (3) 


O 
in which ù is a real number denoting the half-amplitude of 
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the incident wave. All phase angles are referred to the 


wave crest with positive values indicating lag. 


In Eg. 2, m, denotes the dimensionless body mass tensor 
1) 


and is defined by 
= = u 3 
m,=m,,=m,,= m/a (4a) 
uS scu Za "3209 a) 


in which m denotes the mass of the body, £ denotes the fluid 
density and a denotes the characteristic length scale of the 


body. Sisilarly, 


= DONE - - 5 ü 
ma 7 L,, [PAE V HR prf Pat, Tee a (us 


Wyss uuo p Pat, Ryo 17709 EL ¿IRA 
(4e) 
= = - a 
Mr Maem ER a 
where the mcıents of inertia are defined typically as 
= - y! (5a) 
La a 
m 
= - lf (3b) 
Li “Las J yz dn 


in which the integration extends over the complete mass of 
the body denoted by m. 


The excitation force and moment coefficients are defined 
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as 


TT 
- F:(m : 
C; p A Aue) B DO 2 (6) 
Ga 
; ib; . 7 
C. = Ẹmaxj)e 3 T (7) 
Š Jar DE a 
where 6; denotes the phase shift of the force with respect to 
the crest cf the incident wave. A positive value of b; 


denotes a lag. The term F (max), (j=1, 2, ...6), denotes 
J 
the maximum value of the force (j=1, 2, 3) and moment (j=4, 
By 6) acting cn the body. 
In the case of the forces (and moments) due to the 


motion of the body, it is more common to describe the 


hydrodynamic force in terms of the added mass tensor, M . 


1) 
and damping tensor N, . These parameters are defined as 
Fix max) ($ 
Eu. IS 5s s is (8) 
E F: (m ises . : 9 
7 UN;: = AU e 2 £24,556; srl, 2.1.6 ns 
2 uy 
£0*a" X; 
where F (max) denotes the amplitude of the force or mroment 
1) 


Mathe i-th direction resulting from the motion of the 


flcating Ecdy in the j-th degree of freedon. 


The fcrces đefined in Eqs. 6-9 are given in terms of 


surface integrals cf the hydrodynamic pressure as 


ee (t) = y Pj h; ds ) ES 12. .6 d 
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Fi - ff Por h; ds auo e au 


where F denctes the i-th component of wave excitaticn force 


i 
or moment and FP, denotes the i-th component of force or 
mcment arising from the j-th component of body motion. The 


. functions h, are defined as 
i 


D (12a) 
h,= n, (12b) 
h, = nz (12c) 
hy = (d ey) ngz-zny (12d) 
ea n US! 
ng = Xny—(dty) ny Da) 


in which n, n and n denote the three components of the 
X y Z 


outward urit normal vector on the surface of the body. 


Dhexparameters P ,9(9-21, 2, ..-.6), denote the pressure 
J 
on the immersed surface associated with motion in the six 


degrees of freedom and MOT denotes the pressure associated 
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with the incident wave (0) and scattering (7) cf the 


Emeident wave. Finally, in Eq. 2 the dimensionless 


coefficients K denote the hydrostatic restoration fcrce 
i 


(or moment) tensor and are defined as 


nciht _ 2 
K,,” [| ny ds = Me (13a) 


| 
TE NS a3 {{ zn, ds (135) 
S 
| +; ( v'z' Ee )d (13c) 
t47 at e ez eS 
A 
LY. -— i (13d) 
Kye Koy: a EZ ny as 
Koc E us X n, ds (13e) 
S 
| 
Kog = 75 | (xy n, d p) ds (138) 
S 


where A denctes the waterline area. 
W 


Mooring line reactions depend, in general, on tha 
mocring line configuration, weight, shape (catenary), 
hydrodynamic and elastic properties. However, it is 
generally permissible to disregard dynamic effects and 
approximate the reacticns by a linear relationship with the 
six components of body displacements as indicated by Eq. 2. 


The term K' in Eg. 2 represents the dimensionless force 
1) 
(or moment) produced in the (negative) i-direction by a unit 


displacement of the body in the j-direction. Given the 


mooring line configuration and elastic properties, the 
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Spring constant tensor can be evaluated. 
D. REPRESENTATION BY VELOCITY POTENTIALS 


The fluid motion in each of the seven flow situations 
previously described can be represented by a velocity 
potential. Specitication of the velocity potential of a 
flcw completely describes the flow, since fluid velocities 


at any pcint are obtainable through the relation 


q 2,8 7 VO(x,y,z,¢) T 


where qQ(X,y,Z2,t) represents the velocity vector, D(1,Y.,2.t) 

represents the time dependent velocity potential, and 

VG denotes the gradient of the potential. Thus, for the 

Six modes of oscillation of the body in still water, the 

fluid velccity vector is given by 

-iot 
] 


q Oo pz t) =Y y, zt) = Re [V4 21 De (15) 


ko aC 
in which MD denotes the compiex space dependent part 


of the total potential induced by the k-th mode of 
oscillation. The velocity potential associated with tne 
seventh problem, waves interacting with the fixed body, is 
actuaily composed of two components. One component is the 
potential due to the linear incident wave, the space 


dependent part of wnich is given by 


De Las bs k (ht gcns 


S on tosh Kh 


m 





in which 9°=H/2 denotes wave amplitude (half-amplitude) and 
k= 2N/L denotes wave number. 


The seccnd component of the wave interaction potential 
is due to the presence of waves whose characteristics have 
been altered by interacting with the object. This potential 
is referred to as the scattering potential and is denoted by 
the subscript (7). The total potential for the motion cf the 


body in regular waves is expressed through the superposition 


s 
17 
2 Na t) ug 


$ (xy, nc 5 
z,t) E, 


E. THE ECUNLARY-VALUE PROBLEM 


Each cf tbe velocity potentials ( $. k=1,2...7) must 
satisfy the continuity equation, which in the case of 


irrotaticnal flow takes the form 


V* $ (42) = (poo 7 (18) 


In additicn, the pctentials must also satisfy the kinematic 
bottom bcundary condition 


Soe (x,-h, 2) =0 Mo. (19) 


where h denotes the mean water depth, as well as the 
linearized dynamic and kinematic free surface boundary 


conditions exrressed by 
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(20) 
) 


DÈ 2 = 
5,500062 - € $,(X,0,2)-0 Ko v T 
where g denotes the acceleration of gravity. 


While the free surface and bottom boundary conditions 
are exactly the same for all of the potentials, the 
kinematic bcundary condition applied on the surface of the 
object is different for each case. This boundary condition 
results from the assumption that the surface of the object 
is rigid and impermeable; it simply imposes the condition 
that there be no normal component of fluid velocity relative 
to the surface at any foint on the surface, and is expressed 
mathematically as follows: 


For oscillation of the body in the still fluid, 


o SE (21a) 
20s - g2=-ieX,ny (21b) 
- Ta men (21c) 
= = 9, = -ic 8, [(drynz + zn,] (21d) 
E. ds -i& 6- [znx- xnz] (21e) 
E =9,= -i0 &, Uxny - (de) n,] (212) 


where n , n anā n denote the x, y and z components of the 
x Y Z 


unit outward normal vector on the immersed surface which is 
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defined in its mean position by S(x,y,z)=0. 


For wave interaction with the fixed body, 


DY, _ 9 96, - Lq7°k eilkxcosyrkz sin y) 


on ? $n Fes (kh) 


. [n Sinh Ük (hey)] e C (n, cos Y «n, sin z)* 6osh [k(hry)l] (22) 


in which Y denotes the angle between the incident wave 
propagaticn direction and the x'-axis. Equations 21 and 22 
are applied only on the immersed surface as defined by 
S (x,y,2z)^0. 
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Bee GREEN*S FUNCTION SOLUTION 


In using Green's function representation of tne velocity 
potential, sach of the seven unknown potentials is 
considered to be induced by fluid sources which are 
đistributed over the surface or the body in a continuous 
manner. The strength of the source at any location on the 
Body surface is a function of position and is denoted ty 
EN. where (527 95) indicates a point on the immersed 
surface. (Figure 3 illustrates the concept of the continuous 
source distribution.) The velocity potential associated with 
each of the bcundary value problens is then given by the 


surface integral 


$ (1,2) = ff F(§,7,3) G(X,¥,2;4,2, 8) ds (23) 
S cd T 
where G(X,Y,2i $,71,5 ) is referred to as the Green's 
function. its value depends upon the location of the point 
[7,7 ) On the source surface, and the point (x,y,z) in the 
fluid field. The Green's function for the present case, 
which satisfies the bottom and free surface boundary 


conditions, is given by 
G= R + yR! + Gt (24) 
in which 
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Rz2AKG-xY +(9-y (5-2) (25a) 
R' e 4829" + (y rahe a) 46-27 (255) 


The G* term is given by Ref. 3 as 
* co 
G (11235, = 2PY f 
O 


(ure cosh Lu (neh) cosh La y o4)] I, rid ut 


A sith (Ah) - ¥cosk(4h) 


+ in (key) cosh [K(7+4)] cosh [kty+4)] J, (kr) (26) 


k^h -»*h «y 


in which 


r=f(x-§)* +(Z-S)° a 


The tern ur denotes the Bessel function of the first kind, 


cf order zero, and y represents dimensionless wave 


frequency. E.V. indicates the principal value of the 


infinite integral. 


An alternate series form of the complete Green's 


function is given ty 


G= an cosh [K (h«7)] cesh Lk («eJ Y (kr) - (X, 
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oo 


+4), an cos[4, (yrh)]cos [4 ENLA (4, r) (28) 


kzi 
in which rn denotes the Bessel function of the second kind, 
of order zero, K denotes the modified Bessel function of 
the second kind, of order zero, and 44, represents the real 


positive roots of the equation 
Hy tan (4, h) r)»-0 (29) 


In numerical evaluations either form of the Green's 
function may be used with the exception of the case where r 
approaches zero. Here the series form given in Eg. 28 is 
Singular so the integral form given by Eq. 26, which has a 
well-defined 1/R singularity, is used. However, when r is 
not smali the series form is generally used for numerical 
calculaticn kecause it requires less computer time for 


evaluaticn. 


The Green's function representation of the potential 
given by Eq. 24 satisfies the laplace equaticn as well as 
the free «surface and bottom boundary  conditons fcr all 


cases. The scurce strength function peel? 1s unknown 
and is determined by application of the kinematic boundary 


conditions cn the body surface (Eqs. 21 and 22). This 


results in tke following integral equation: 


/ 
ir, f ($,%, $) xc ios = $ Getz) (30) 
| pce T 


in which ?G/3, denotes the deriviative of the Green's 


function in the outward normal direction, which is given by 
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BS (X42, 87,9) 2c uz» s, 7n (30a) 


HE ren represents the outward unit normal vector at (x,y,z) 
and where (x,Y,Z) now represents a control point on 
Sek, Y,2)=0. 


For the seven problems, ail quantities in Eg. 30 are 


mew except for the source strength functions f ($,72,9); 
K 
(k=1,2, ...7). However, the integral equation cannot be 


solved analyticaliy for an arbitrarily shaped body, and a 
numerical approach must be adopted. The numerical approach 
developed herein involves discretizing beth the source 
strength function and the pody surface beundary condition 


functions, Sale) 


See DISCRETIZATION Or IHs PROBLEM 


Ihe technique previously used [Ref. 1] to discretize 
tne integral equations given by Eq. 30 is to divide the 
Esurce distribution into a grid of N quadrilateral panels, 
and to assume that the source strength has a distinct, 
constant value over each panel. Hence, for any of the seven 


ficw problems of interest here, £ denotes the strength of 


J 
the source everywhere on panel j, (j=1,2...N), for the k-th 


problem (k=1,2 ...7). The area of panel j is denoted Ly 
AS; and the induced potential anywhere in the fluid thus 


msultsS freo a summation of the effects of all the Source 
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panels. Figure 4 depicts the discretized source surface and 


resulting potential. 

Additionally, the normal boundary condition 
functions nez) ig uort-E5g. 30 are discretized 
by applying the kinematic body surface boundary condition at 


only M selected points on the surface. In practice, the 
humber of boundary condition ("control") points are chosen 
to be equal (M=N). For convenience, the control points are 
chosen to lie at the centroids of the source panels. Thus, 
for each discretized integral equation, a system of N linear 


eguations in N unknowns results: 


N 
/ | 9G (i, 5,20 £59) 45 - en 207 
— gaur d >= 2 : 
41f 2. n N SL Jk 95) 2 NEN) 
in which (x ,y.,Z.) represents a controi point i on the body 
B TT 


a race, at which the boundary condition is being apgliea, 


and S, denotes the value of CI date pont i. Equation 
i 
31 is expressed in matrix form as 


[x] (d = EN EE 0T (32) 


where [&] is an NxN matrix composed of elements 
5 DE (X.Y; z;.£,2,5) 33 
eai n A) E Pic ( ) 
31725, 
The vectors tt] and ij represent the N unknown source 


Strengths and the kinematic boundary condition applied at 


mimes N ccntrcl peints, respectively, for the k-th flow 
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problem. It should be noted that the matrix [oc] depends 
only on the body configuration, water depth, and period of 


the motion ard is the same for all seven cases. 


Ihe discretized potential which is induced by the N 
Source panels is specified by a matrix equation similar to 


Ec. 20: 


(9, z le | (R) 33 2007 (34) 


The martix [fB] is composed of elements 


N 
an A f G(s,2,5; 1; 1,2) d's E 
J=/ "AS; 


As in the case of matrix [c], matrix [(3] is the same fcr all 


seven flcw cases. 


Integration of the Green's ‘function and its 
derivatives in Eys. 33 and 35 poses no great difficulty. It 
is most convenient to separately integrate the singular 
terms. Special considerations in performing these 


integrations and necessary rormulas are provided in Ref. 2. 


A ccmputer solution of Eq. 32 is accomplisned using 
Matrix inversicn and multiplication, yielding the previously 
unknown values of the source strengths. Values of velocity 
potential at the control points (panel centroids) are then 
obtained directiy using Eqs. 34 and 35. Fluid velocities at 
the contrcl points are similarly obtained after first 
differentiating the elements of the [Q] matrix with respect 
to the x,y and z directions. Dynamic pressures are 
calculated frcm the linearized form of Bernoulli's equation. 
Finally, the forces acting on the object are determined by 
pressure integraticn (cr more precisely, summation) over the 
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surface area of the object as indicated by Eqs. 10 and 11. 


2.) Method of Linearly Varying Source Strength Element 


An alternative method of discretizing the Green's 
function integral equations, proposed by Webster [Ref. 10] 
for unifcrm unbounded flow past a submerged body, utilizes 
triangular scurce panels over which the source strength 
function is considered to vary linearly. In this way, the 
discontinuities or "jumps" in source strengths at adjoining 
Panel edges (which existed in the previous scheme) are 
eliminated. Also, since source strengths are not held 
constant over the extent of each panel, the discretization 
results in a less crude approximation of the continuous 
strength function. An additional "smoothing" effect on the 
velocities induced near the body is achieved by submerging 
the source panel surface inward a selected distance from the 
actual body surface, although the kinematic body boundary 
condition is still applied at points on the true body 
surface. It is noted that while some smoothing effect can 
be accomplished Ey submerging the sources, this is not a 
necessity as long as the Singularity is accounted for 
properly. 


Figure 5 depicts the linearly varying source 
strength function for a single triangular panel. The local 
coordinate system for each panel is defined such that two 
triangle corners (denoted by a and d} lie on the local 
X-axis, and the third corner (b) lies on the positive 
y-axis. Appendix A gives the transformation between the 
body coordinate system and the local coordinate system for a 


given parel. 


In panel coordinates, the linearly varying source 
strength function may be written 
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L(8,7)= +85 +37 (36) 


in which &,f and Y denote constant coefficients and the 
tildes indicate lccal coordinate values as defined in Fig. 
5. The values of the source strength at the three ccrners 
are 


m = F(a,o) (37a) 
Fo 
f 


- (o, t) | (375) 


a) (370) 


T 


Sigstitumiyon Oc has. 37 into Eq. 36 allows the 
constant coefficients to be expressed in terms cf the 
triangle geometry and the values of the strength at the 
SeEners: 


e =(f,-E)/(a-d) e 
« z (afi -4 £,)/(a-4) (38b) 
3 - [(a-4) f, -a £j rdf] Joa) (38c) 


Baus, the  ccntribution to the potential at (x,y,z) due to 


the single panel of area AS is denoted by AQ(x,Y,2) and is 
given by: 
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Www / Ba A 
apaa f (e «8 $17) e (Sr x zz) ds (39) 
AS 

where $ and 7| denote local panel coordinates and (x,y,2) 
denotes a general point expressed in terms of Local 
coordinates. 

Substituting Ego. 39% into Eq. 39 and collecting 
terms allows the contribution to the petential to be 


expressed in terms of the values of the source strength at 


the corners cf the panel: 
Ab(K,7,2) =f AG, +f, 46, + Fy Ad, (40) 
n which 
A>,= en ds dig 
A®,= f. t (£7 .X, 7, Z) ds m 
Ady = “ff % ee CAG, 7; Lizads (41c) 


For convenience Ad, od, and AA will be referred to as 
Simply the “corner integrals" associated with corners a,b 


and d, respectively, for a given trianqular panel. 


The velocity induced by a single source panel is 
ottained Ey taking the gradient of Eq. 40: 


AQ Cie) = f, V(Ad,)+ fF V(ad,) + Fy VA dy) (42) 
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The gradients of the corner integrals are given by 


_ (( - bd tb Fede a 
vag) = f. um b@-d) VARTA MEO 


V (Ad) 7 f. I VG($,2;:3,5,2) ds (43b) 


a een a T 
V (adj) i ).. 4 fb(a-4) Vet 452)45 i 


The evaluation of the three corner integrals (Eq. 
41) and their gradients (Eq. 43) must be accomplished in 
local ccordinates for each panel. In performing the 
integrations, it is most efficient to separate the Green's 
, function into its singular and non-Singular components, and 


BEC OO une ee conoce: 
sh, f va (rids +f eds +f 4 (GA)ds Lm 
Ah, = ff m (Vds ef wa (Ves ff Mb (eds wen 
bd, ={) Wa (Rds + f( va («Das eff. ES (440) 
V(Ad,) = f NL V(/8)4s « oe) as : (were nac 
Tabi S, Werd, m URS Mer” 
v(58)) - )| wav(V0ds » (f wu v CA gdse favie)ds S 
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Sur 


in which i , Er and à. are linear functions in $ and f aná 
a 


are given by 


a T bd +b F +d Y 
Mae an bla-d) nn 


M ko SS | (46b) 


Wy = ba-b$-a7 (46c) 
4 T b (a -d) 


Àppendix B provides foraulas developed by Webster 
[Ref. 10] for the closed form evaluation of those integrals 
in Eqs. 44 and 45 which involve derivatives of 1/R and 1/R! 
as well as similar formulas developed by Yeung and Bai [{ Ref. 
11] for integrals of 1/R and 1/R'. The integrands are 
Singular, and it is therefore necessary to evaluate the 
integrals analytically. 


Integration of the G* components of Eqs. 44 and 45 
is made possible by assuming that the G* function and its 
gradient vary linearly over the triangular panel, (as was 
done for the source strength function). This is a 
reasonable approximation since G* is a well-behaved function 
that varies slowly with position on the body surface. The 
wave length of its oscillation is approximately the same as 
the wave length of the incident wave, and this tends to be 
large relative to the panel dimensions in typical 
applicaticns. 


With G* and WG* expressed as linear functicns of 


their values at the panel corners, the terms in Egs. 44 and 


45 involving these two functions may be expressed as 
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f Ma Gds = n (a,b, d, GG? EL) (47a) 
ff, Wee" ds = fy(a,bd, &* cf, c7) (47b) 
Í, wjG'ds * I4 (a, bd, CG. an) (47c) 
$ Wa VG'ds = T (a,b d ; EON VG, VG) (48a) 
Da Wa VG "ds = í, (a, bd, a. uve) (485) 
I, wi Ve*ds = L (a, b d, IET IGF, 764%) (48c) 


daere ,) 1 , 17, D" 1", It, represent the integrals cf the 
a b d a b d 


product cf the two linear functions over the triangular 
area. The integrals in Eqs. 47 and 48 were evaluated using 
the formulas given in appendix C for integrating the product 
of two linear functions over a triangular area and the 


resulting integrals are given in Appendix D 


The total potential and velocity at any point in tae 
fluid is a result of the combined effects of all the source 
panels. For a grid composed of M panels, the potential is 
given by 


mM 
$0,527 X. la 2d, (01) + 09, (64,0) e fL Ad, (ott) (49) 


and the velccity vector is given by 


E M | 
q (x 14)Z) S 2 ha M (^ ak 12) + f DÍA Ao 2)) B fa ve $ NY y 2) (50) 
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in which f y E e  , denote the strengths at the three 
a D d 
K K K 
corners of the k-th panel. 


It should be emphasized that the corner integrals 
and their gradients are evaluated in their own individual 
panel coordinate systems. Transformation back to the global 
cocrdinates, described in Appendix A, must be accomplished 


pricr to performing the summations of Eqs. 49 and 50. 


The summations given in Eqs. 49 and 50 can be 
re-arranged by reversing the order of the Summation. That 
is, if coefficients or the strength at a given node are 


first summed, then tne summation may pe written 


N 

(x,y,z) = 2 F AQ: (X,4,2) (51) 
J= 

_ N 

q(x4y2)= 2 f; V(A0;(xy,D) (52) 


| 
nach £  denotes tae source strength at corner node j 
J 
e-1,2...N and ^d. represents the sum of the corner 
J 


integrals associated with node j, for all panels which share 


ncde j as a mutual corner. A graphic representation of the 


meaning of DO. eoeprotrded 1 Fig. 6 for Clarity, since it 
J 


is difficult to exprəss mathematically. As inđicated in 


- 


hag. 6 Ad (X,Y,2) denotes the contribution to the 
J 
potential at (x,y,z) fron the panels surrounding node j. 


Each of the surrounding panels may be considered to have 
unit strength at the corner coincident with node j anā zero 


at the otker two corners as depicted in Fig. 6. 
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The kinematic boundary condition at the body surface 
may now ke applied at N selected points. It is convenient 
to use points on the surface which lie directly outward from 
the corner nodes of the submerged source surface. In 
practice, E is best to initially define the nodes cn the 
actual body surface, and then establish the points on the 
submerged grid at scme pre-selected submergence distance 


alcng the corresponding inward normal vectors. 


Application of the kinematic boundary conditions for 
each of the seven prcblems of interest results in the 


discretized versicn of the integral equation: 


Me ce» 


EU - 9. 


where the repeated index j indicates summation. The term 
denctes the sum of the corner integrals for node j, 
evaluated at contrcl point i. The derivative with respect 


to the outward normal vector is obtained by 


AQ; i Aan $ 
E ad = V (AQ; ;) KE (54) 


where n, is the unit outward normal vector at the ccntrol 
T 


point i, cn the body surface. 


Equation 53 may be expressed more simply as 


[=] 16) = 49. ERU T i) 


where the NxN matrix, [&], is composed of elements 
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E 9 (20,,) (56) 


Equation 55 is solved by computer, (forn ksl eaa l) 
using matrix inversicn and multiplication, yielding the 
unknown source strengths. Substitution of the source 


Serengths into Eqs. 51 and 52 with evaluation at (X.,yY.,2z)) 
AL 

(i=1,2,...N), gives the values of the potential and the 
fluid velocities at the control points on the body surface. 

Tke dynamic fluid pressure is determined using the 

linearized form of Bernoulli's equation. For oscillation of 


the body in its six degrees of freedom, the pressures at the 


surface ncdes (denoted by i) are given by: 


E r. -igt 7 
PL. = Relivo >, © ] en... (57) 
; gt 
nen: =Re LiPo (do; ty. Je ] [29) 
i t t 
where denotes fluid density, and P and $. denote 
i i 


: \ 
dynamic ¡pressure and potential, respectively, at node i 


(denoted by the subscript i) due to motion in thea k-th mode 


of oscillaticn. 


The pressure integrals given by Egs. 10 and 11 
cannot te utilized directly to obtain total forces and 
moments, because only discrete values of pressures are 
available. However, it is consistent with the numerical 
Scheme tc assume that the pressure varies linearly over each 
triangular panel. Based on this assumption, the pressure at 
any point on a surface panel is expressed as a linear 


function of the panel geometry and the values at the corner 
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nedes, which are given by Eqs. 57 and 58. 


The forces and moments are then calculated using 


mM 
la) = f Fk h; ds 3 Dir een SA 


ot “AS 


ei 
7 2, if D ds eG = 
|^! AS: 


where F, denctes the i-th component of wave excitaticn force 
1 
or moment and 4 denotes the i-th component Of LOECE OF 
1 
mcment arising from the k-th component of bcdy moticn. The 


summations cf Eqs. 59 and 60 are carried out over the M 
panels on the actual body surface. The functions h, are 
i 


given by Eg. 12 and are functions of position on the 


immersed surface. It appears to be consistent with the 
numerical scheme to approximate these functions as linear 
functions so that they may be defined at interior points by 
their values at the corners of the triangles. Thus the 
Macegrands of Eqs. 59 and 60 contain the products of two 
linear functions. The integrals are evaluated in a 
straight-forward manner using the integraticn formulas given 
an Appendix C . 
The resulting forces F, Gk 1726.20) “Sand F 
ik 1: 
E21, 2...6) are substituted into Eqs. 6, 7, 8 and 9 to 
evaluate the excitation force coefficients C (j=1,2,.-..6), 
and the added mass and damping coefficients, M , and N, , 
i i 


(i,j=1,2,...6). Finally, the equations of motion for the six 


degrees cf freedom (Eq. 2) are applied to determine the 
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ccmplete dynamic response of the body. 
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A computer code based on the numerical method outlined 
was developed to calculate the hydrodynamic coefficients as 
well as the dynamic response of a floating body in waves. 
Using this program, example calculations were made for two 
simple gecmetric shapes: a vertical circular cylinder which 
extends from the bottom and passes through the free surface, 
and a semi-inmersed sphere. Computer produced drawings of 


the grid configurations used are shown in Figs. 7-10. 


The vertical circular cylinder represents the only 
three-dimensicnal geometry for which a closed-form solution 
exists (of the type of interest in the present study), and 
therefore, it is of interest for making a comparison of 
numerical results computed by the present method. MacCamy 
and Fuchs [ Rei. 12] have developed a closed forn solution 
er the herizontal force acting on the cylinder and this 
BESUIt is plotted in in Fig. 11 in the form of the 
dimensionless amplitude of the force versus the wave length 


parameter, 2Mfa/L. 


Corresponding results computed by the triangular panel 
method using grids of two different finenesses are shown on 
Erg. 11 for comparison with the results of MacCamy and 
Fuchs. These results indicate rather rapid convergence of 
the solution, particularily at small values of 2Ma/L. The 
Slower ccnvergence (or greater error) at the larger values 
of 2MWayL apparently indicates the inaccuracies associated 
with the assumption of linearly varying source strength. As 
2MWa/L Leccnes large (the wave length becomes small) the 
Variation of the G* part of the Green's function with 
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distance beccmes more rapid. Thus, as 2Ta/L increases, a 
point is reached where the assumption of the linear 
variation of the Green's function becomes appreciably 
Envalid. Of course, as the fineness of the grid is 
increased the numerical results tend to agree witn the 


clcsed-fcerm exact results up to higher values of 2Ma/L. 


Havelock [Ref. 13] has presented results for the added 
mass and damping coefficient in heave for a semi-immeéersed 
sphere and these results are presented in Fig. 12 fOr 
comparison with numerical results computed by the present 
method. ihe figure shows the heave added mass and damping 


coefficients, es and ER respectively, EOD the 


semi-immersed Sphere along with the numerical results 


computed using grids of two different finenesses. These 
results indicate that the numerical results are tending to 
@Genverge althcugh even for the finer grid the results are 


not yet converged. 


Moreover, given the damping coefficient, the wave force 
coefficient can be computed by use of the Haskind's 
relations which, in general, relate the force coefficients 
Ew the far-field solutions of the corresponding radiation 
problem. For the case of the heaving motion of an 
axi-symmetric body a closed form relationship exists between 
the heave damping coefficient and heave excitation force. 
Thus, Havelcck's results for the damping coefficient in 
heave were used to compute the heave excitation force and 
Ense result is presented in Fig. 13 for comparison with the 


numerical results. 


The numerical results presented for comparison with 
mawelock's results in Fig. 13 show a trend toward 
convergence and, in fact, the results corresonding to the 


finer grid indicate adequate agreement. 
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Finally, the dynamic response in heave and surge for the 
case of the semi-immersed sphere are shown in Fig. 14. 
These results, showing the dimensionless response in terms 
of the respcnse amplitude to wave amplitude ratio, indicate 


very little difference between the the two grid finenesses. 
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À new nurerical procedure based on the use of triangular 
source panels has been developed for computing the 
interaction of fixed or floating bodies with waves. 
Numerical results based on this procedure ccmpare well with 
existing results, and, in general, the procedure appears to 


converge to the correct solution. 


The value cf the triangular panel procedure as compared 
to the procedure based on quadrilateral panels of uniforn 
strength remains uncertain. This will be determined cnly 


through experience in application of the two methods. 
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APPENDIX A 


COORDINATE TRANSFORMATION 


Putas following, the coordinate  transformaticn is 
developed which relates the local coordinates attached to a 
particular panel, as indicated in figure la, to the global 
coordinates with origin located at the mean water level. The 
global system is denoted by 0(x,y,2) and th local 


Coordinate system is denoted by 0(X,Y,2)- 





Z 


Figure 1a-Coordinate System Locations 


The location of the origin of the local coordinate 


System expressed in global coordinates is given by 


pz Fe (X7 X1). (ara (Ya), (ZA Fc (z.-zi]] ' 
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in which (x 203. X Z and (x Z are the global 
( EE d ( d 2 ( Ae ay a’ g 


System cocrdinates of the triangle corners a,b and d, 


respectively, and where 


+(Za-Zi)(Z6-24) 
(Xa- X4)* ed) (Za u) 





The unit vectors in the x,y and Z directions of the 
Pecal ccordinate system are given in terms of their 


respective global system components as 


ü 


Uy i P Ux, d ty, k 


el 


= Uy + Uy J + 4y K 


I 


Uz i + Uz j + Uz, k 


Mu Gay —-y 


, j and K are the unit vectors in the X, Y and Z 


y 
i 


in which 
@ucections, and i, j and k are the unit vectors in the x, y 


and z directions. 


The terms Ue, UW, ...4 are defined as fcllows: 


X, x, 2j 


Ux = Q7 X/ \ dal 
Ux, = (a7 f2)/ 1 dal 
a 
Uy, = (4 X0)/1 ob! 
Uy, = (Y, -1.)/10bl 
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ae za obl 
u x Uy, Uy Ux, Uy, 


"e = Ux, Uy, Ux, Uy, 


| dal XA Y 44)" HRS ZA] y 
lob| =. /&- X? Hu 1° He zu” 


The transformation matrix is then defined as 


where 


and 


Ux Ux dx- 


T] > Uy Usa Uy | 


Transformation between global coordinates and local 


coordinates is accomplished in the following manner: 


X E Á- X, X T X X. 
O Y'zl[T] 17 0, 
Z ZZ; AE > Ze 
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T 
in which [T] denotes the transpose of the matrix [T). 
Similarly, a vector quantity, F, evaluated in tne local 


coordinate systen is transformed back to global coordinates 
using 


= [1] 


ON 
An 


me waoich E , F and F are the x, y and Z components in 

xX y Z 

local coordinates and F, F and F are the correspcnding 
X y Z 


components in the global systen. 
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APPENDIX B 


FORMULAS FOR EVALUATING INTEGRALS INVOLVING 1/R AND 1/2! 


The terms ot Eg. 44 which involve 1/R and 1/R' are 


evaluated as follows, using the panel geometry defined by 
Big. Tb. 


y,” 





: x.8 
(d, c) (2,0) z 
Figure 1b- Panel Geometry 
-bd +b -bdsb$ cd | 
er 4 b(a-d)R — 2o (a-d)«n Cen TRAR 


(brad) Eu /c? - Ey * a (e rad)(2? - Zn 


+e(22- Z') ul - (25»'-d(?«2)] x 
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ba-b §- any J Dr 
f Ambla-dR — 2b(a-d)un {-292 -a Rad 


+(b*+ ad) Fy, /e? - Fp * [222'- (brad y] Ty, /e 
eZ -£,)1 E A z Sey ¢ 

f reeds = a f eR ca ls, [ c* -d Ep /e? 
+ Lapp a (Z^-A)]EÀ /c — CE - 2) 

- Eag»'- d (2? +22) rh, Je 


in which the following definitions are applicable 


n= (X-a) + "n rz* Caw ee 
2I ES Cice E 
= X- + (Y-b) E area -b 


= Ud ++ 3? ez” 


P= bX tağ -ab P odi - 

Paz aX -b7 -a* A! dX - bs of’ 

we = ax 2.09 Beer 
^ i o aa 1a = (X-a) 2 
== 127 a ET UT a] 


+ tan E8523] — tun Led? p UM 
o Fa 


nz eae nas 1 (X-d)] 
as 7 log (Crst A) (en & Ay] 
Tj, = leg[(era £A4)/(ev, t7] 
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Rad = Va (K-2) rg ($-d) Moe 


-— . 7 A 
Fab E Fy" = f Lab ” ae) 
E / _ S — 
Fb A E a pd a 


Formulas for evaluating the teras of Eq. 45 which 
involve tke gradients of the 1/R and 1/R' integrals are 


given by Webster [ Ref. 10] as 


LAA (de = Í7LI+ETP_ dy 
SGT alleles = (LEE Bon 


e 
2312, - 225] 


p 
+ Lean edad) 


b -b$-a75 _ EE ERS m. 
|. ED BUR) ds = n EL * B n 
+ 
* L^: bo 227 rf, — Plat $ 
LED B 
J), Gm Rds = ER 7 4 ir 4) - Sr z,- 211] 


(| -bdtb$td] ZUR) ds = - $2 tz Lad we a Es 





a ue ew) or ba 
A A d ate, + _ AP Thy 
E +L Ee er 2o be (a- A) 
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I. 


BE d. - (oa) 


Lo E 
- 3/85 - UU Men. PO 


ba-b$-4 
e de - 12 he Se 
rad [E ande a) 
" A x = (Hees Sle —Y7) Sla- 15 
A / 
+ [2£ n IRE - [E 1 roe Ine} 
AvR)ds = x mC LIT MBEIII Fe 
-d2 Da4 $ 


» 


j 


meak kerdh 5 


bla-d) | |-P7- Gn ay 


ba=b5-47 2, 
o a-b? -4 a Rds -=r 
i +CZ Lo, eur 
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IPS TX 


FORMULAS FOR INTEGRATION OP LINZAÀR FUNCTIONS OVER TRIANGULAR 
AREAS 


All Lormülas given are applicable for a triangle 


positioned as shown in Fig. 1c, with a>d and b>0. 


(0,6) 


> 


o (d,o) (a,0) 


Figure 1c- Triangle Gecmetry 


If £(x,y) and g(x,y) are any linear functions expressed 


as 
F(X) = +Bxr37 
~ i Fo p 
INN= xX FBX Y 

and the function values at the triangle corners are denoted 


50 





by zu E fa Ber I” e then: 





= afa-dh ez aq, - q.d 
a-d | 

B= bate g'z dc d 
a-d 

Y «(a-d fs-at, rf 3-6-4)13,- 241 eda, 
b(a-4) b (a4) 


Ende the integrals of £(x,y)ds and f(x,y)g(x,y)ds are given 
by 


f. FR Nds = Vs ar fa) w 


fj £(x,7) 9(X, 7) ds = aa Wir (Ba roca )w; 


+(X ta ZYW (XB+ BIY Wa +BB We FY! We 


where 
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we pe ae 
Du =f, "PB = Hier - 


Ws = j X ds = (a?+d rad Ya-4) b 
de 


W; = N y ?4 s : "^ 





APPENDIX D 


INTEGRATICN FORMULAS FOR INTEGRALS INVOLVING G* AND VG* 


oo e mnt onor a linear variation of G*(£775xX,Y,2) 
and the subsequent 
formula for the integral of the 


of two linear functions, 


over a given triangular panel, 
application of the product 


given in Appendix C, leads to the 
following expressions for Eqs. 47 and 48 which involve the 
G* term. Geo erm onig: ID foc applicable triangle 
gecmetry.) 


( —bdrbird% RE: 
), regn) 08 mm Aare 


AS 





END In pe eer 
)\. aa Gee G ds 4T (a-4)* Ce r C, B, FC, 8, | 


> 


a * "D | 
f. 4 Tb Es 471 b(a-i) LC, D, +C, Pa C C3D,] 





E tonE 
AN Sy 


Es [(a-d)G*-a G} ra GT] /> 


C,= Lac -dG*] 
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- b(2a-d NES E 
M = Bead) = & A, = b, 


AT R 2% 

B, =(atd)b Bas toe B, = ap 
24 s 
bOd ta) p, = bt p 
24 24 


The integrations involving the gradients of G* similarly 


result in the following expressions: 


A 


(| ~bd +b d3 a | 
NES etd = ebenen.) 


batbs-a7 
fj mL 4b (ad) Zar na rH (a-a)* - Kec) 8, (V 62) 8, «Vc, ] 


| 4i b VG'ds - DURO dl) [(vc))5, AE DEN DA 
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Fig.-1- Problem Definition 
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Fig. 2- Definition of Body Motion 


57 





Fig. 
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Fig. 4- Uniform Source Strength Discretization 
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. 5- Panel of Linearly Varying Source 
Strength 
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Figure Jf OURSTER@Or A 13 "NODE HEMISPHERE 
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Figure 8 - QUARTER OF A 25 NODE HEMISPHERE 
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Figure oe OUR eR OF A 49 NODE CIRCULAR CYLINDER 
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Luc se CUART ER OF AN 30 NODE CIRCULAR CYLINDER 
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Added Mass and 
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Fig. 13- Heave Excitation Force for a Floating 
Hemisphere 
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